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Abstract 

A  new  class  of  potential  suitable  for  modeling  the  adsorption  of 
water  on  different  metal  sites  is  described.  The  new  potentials  are 
simple  in  form  and  convenient  for  use  in  computer  simulations.  In 
their  real  space  form  they  comprise  three  parts:  a  pair-wise  sum  of 
spatial  anisotropic  12-6  potentials,  a  pair-wise  sum  of  isotropic  short 
range  potentials,  and  an  image  potential.  Two  modifications  of  the 
potential  are  developed.  In  the  first,  the  anisotropic  potential  acts 
only  on  the  oxygen  atom  and  not  on  the  protons.  In  the  second,  the 
potential  acts  on  all  the  atoms  of  the  water  molecule.  In  practical 
calculations  it  is  convenient  to  transform  the  potential  to  a  reciprocal 
space  form  in  the  manner  described  by  Steele  [Surf.  Sci.  36,  317 
(1973)].  Adsorption  of  water  at  top,  bridge,  and  hollow  sites  on  (100), 
(110),  and  (111)  surfaces  of  Pt,  Ni,  Cu,  and  A1  were  studied  using  two 
fitting  parameters  and  the  results  compared  with  previous  theoretical 
calculations. 


1  Introduction 

Understanding  the  properties  and  dynamic  behavior  of  water  in  the  vicin¬ 
ity  of  metallic  surfaces  is  of  fundamental  importance  in  electrochemistry, 
catalysis,  and  the  study  of  corrosion.  A  variety  of  experimental  techniques 
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exist  capable  of  providing  data  about  surface  water,  for  example  those  de¬ 
scribed  in  the  review  of  Thiel  and  Madey  [1],  surface  infrared  spectroscopy 
by  O’Grady  [2]  and  Melendres  (3],  and  the  X-ray  technique  described  by 
Toney  and  coworkers  [4].  Computer  simulations  play  an  important  role  by 
providing  detailed  microscopic  insight  about  the  water  multilayers,  which 
is  currently  unavailable  from  laboratory  experiments.  A  matter  of  primary 
importance  for  performing  a  computer  simulation  is,  of  course,  the  knowl¬ 
edge  of  the  relevant  interactions.  Our  purpose  in  this  paper  is  to  describe 
a  new  simple  potential  function  useful  for  molecular  dynamics  simulations 
of  metal-aqueous  interfaces.  The  weight  of  experimental  and  theoretical  ev¬ 
idence  to  be  surveyed  briefly  later,  supports  the  adsorption  of  water  with 
oxygen  atom  down  on  the  top  site  of  a  metal  surface.  In  contrast  a  pair¬ 
wise  sum  of  12-6  Lennard- Jones  (LJ)  potentials  predicts  the  binding  energy 
order  for  adsorption  as  hollow  >  bridge  >  top  sites.  It  is  also  clear  that 
the  electronic  properties  of  the  metal  play  a  critical  role  in  the  interaction 
with  water  which  is  weakly  chemisorbed  rather  than  physisorbed,  and  that 
the  chemisorption  bond  directs  adsorption  to  the  top  site.  The  key  idea  in 
our  development  is  the  introduction  of  a  simple  angular  dependence  into  the 
traditional  12-6  Lennard-Jones  potential  to  mimic  the  directionality  of  the 
chemisorptive  bond. 

In  the  earliest  studies  the  water-surface  potential  was  assumed  to  have  two 
components,  image  terms  to  describe  the  polarization  interaction  of  molec¬ 
ular  charge  with  the  conduction  electrons  and  12-6  Lennard-Jones  terms  to 
account  for  Pauli  repulsion  and  the  dispersion  attraction  interaction  with 
all  electrons  in  the  core  and  conduction  band.  We  will  refer  to  this  latter 
part  as  the  core  term.  Earlier  efforts  [5,  6,  7,  8]  described  the  interaction 
of  a  molecule  and  the  conduction  electrons  of  a  metal  often  by  truncating 
the  image  potential  to  the  nearest  term.  In  practical  computations,  the  core 
potential  part  often  appears  in  its  integrated  form  which,  within  the  first 
order  approximation,  is  a  10-4  or  9-3  potential.  These  models,  though  intu¬ 
itive,  fail  to  predict  orientational  properties  of  water  molecules  chemisorbed 
on  the  surface  of  a  transition  metal.  Existing  experimental  data  [1]  for  the 
work  function  change  suggest  that  the  water  dipoles  tend  to  orient  perpen¬ 
dicular  to  the  surface  rather  than  parallel  to  it.  In  the  case  of  low  index  faces 
where  the  metal  atoms  are  explicitly  taken  into  account,  the  above  mentioned 
water-wall  potential  leads  to  a  conclusion  that  the  preferred  adsorption  is  on 
the  hollow  sites  of  the  surface  [5].  This  contradicts  the  weight  of  evidence 
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that  an  adsorption  on  top  of  the  metal  atom  is  the  preferred  site,  as  indi¬ 
cated  by  both  the  experimental  observations  [9]  and  the  quantum  mechanical 
calculations  [10,  11,  12,  13,  14].  It  is  therefore  a  challenge  to  find  a  water- 
metal  surface  potential  function  that  is  simple  in  form  and  still  describes  the 
observed  preferences. 

Extended  Huckel  molecular  orbital  calculations  of  a  water  molecule  on  top 
of  a  cluster  of  five  platinum  atoms  [10]  were  used  by  Spohr  and  Heinzinger 
[15]  and  Spohr  [16]  to  fit  the  water-surface  potential  to  a  set  of  exponen¬ 
tial  functions.  Berkowitz  and  co-workers  [17,  18]  have  prescribed  poten¬ 
tial  functions  for  water-Pt(lOO)  and  water-Pt(lll)  interfaces  based  on  the 
Heinzinger-Spohr  form,  but  incorporating  the  lattice  symmetry  and  corruga¬ 
tion  explicitly.  These  potentials  usefully  describe  the  chemisorption  bonding, 
but  not  the  image  contribution  or  it’s  contribution  to  long  range  interfacial 
electric  fields.  In  a  quite  separate  approach  Siepmann  and  Sprik  [19]  have 
modified  the  bond-angle-dependent  potentials  developed  for  covalent  solids 
[20]  by  the  use  of  overlapping  Gaussian  charge  distributions  centered  at  the 
metal  atoms,  as  a  compromise  between  the  image  charge  formalism  and  a 
microscopic  description  of  the  electronic  structure.  These  more  complicated 
potential  models  are  in  qualitative  agreement  with  more  elaborate  ab  initio 
calculations  and  existing  experimental  data. 

2  Potential  Energy  Function 

The  physisorption  of  molecules  on  insulators  has  been  usefully  modeled  in 
many  cases  using  a  pair-wise  sum  of  atom-atom  Lennard-Jones  12-6  poten¬ 
tials  and  r-1  Coulomb  potentials.  The  extension  of  this  model  to  a  metal 
would  describe  the  molecule-core  electron  interaction  with  a  Lennard-Jones 
(LJ)  potential  and  the  molecule-conduction  electron  interaction  with  an  im¬ 
age  potential.  In  practice  the  Lennard-Jones  potential  also  includes  Pauli 
repulsion  and  attractive  dispersion  interactions  for  all  the  electrons  making 
up  the  metal.  This  does  not  work  for  water,  or  for  that  matter  the  adsorp¬ 
tion  of  Xe  atoms  on  Pt  [21].  Many-body  effects  other  than  those  giving  the 
asymptotic  image  interaction  are  expected  to  be  important  for  the  adsorption 
on  metals.  It  is  not  correct  to  treat  metal  surface  as  a  simple  assembly  of 
uncorrelated  atoms,  nonadditive  interactions  should  be  incorporated.  How¬ 
ever,  to  obtain  a  concise  potential  in  a  form  which  is  convenient  for  adjusting 
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parameters  to  be  useful  in  computer  simulations,  it  would  be  'very  helpful  if 
one  could  treat  the  many-body  contribution  in  an  approximate  but  effective 
manner.  In  this  section,  we  show  that  this  possibility  indeed  exists.  Intro¬ 
ducing  angular  dependence  into  the  traditional  LJ  potential  can  make  the 
top  site  most  favored  for  a  water  molecule  in  an  oxygen  down  geometry. 

The  two  potentials  described  below  are  labeled  A1  and  A2.  The  later 
being  the  all  atom  model.  In  both  there  is  an  LJ  potential  modified  by  incor¬ 
poration  of  angular  dependence.  Instead  of  spherical  symmetry,  the  repulsive 
and  attractive  parts  have  the  symmetry  of  an  ellipsoid  of  revolution.  Several 
interactions  contribute  to  angular  dependent  terms:  first  is  the  chemisorp- 
tive  bond  interaction  involving  mixing  of  orbitals  from  the  metal  and  oxygen, 
another  comes  from  the  interaction  between  the  permanent  dipole  moment 
of  the  adsorbed  molecule  and  the  induced  moment  of  the  solid  atom  [22]. 
Repulsions  due  to  overlap  of  electron  clouds  at  short  separations  are  also  con¬ 
sidered  in  this  modification.  In  addition,  extra  terms  proportional  to  r_n  sure 
introduced  to  mimic  short-range  bonding  that  pulls  the  oxygen  atom  closer 
to  the  top  atom.  These  angular  dependent,  short-range  interactions  allow 
better  modeling  of  the  adsorption  of  water  on  metal  surfaces.  In  simulations 
where  metal  atom  dynamics  and  surface  reconstruction  are  not  considered,  it 
is  useful  to  explicitly  incorporate  surface  lattice  translational  symmetry  us¬ 
ing  the  method  first  proposed  by  Hove  and  Krumhansl  [23]  and  subsequently 
universally  adopted  [24,  25]. 

In  both  models,  the  potential  for  one  molecule  interacting  with  the  metal 
is  written 

Vw—met  ~  Vw—ccnd  "t"  core  (1) 

where  the  first  term  represents  the  water  molecule-conduction  electron  po¬ 
tential  and  the  second  term  the  molecule-core  electron  contribution.  As 
usual,  we  approximate  the  first  term  by  a  classical  image  potential.  The 
image  charge  qim^e  of  a  classical  point  charge  q  is  located  at  the  symmetrical 
position  below  the  image  plane  and  has  the  magnitude  [26] 

?iraa«e  =  .  7  (2) 

I  t 


The  ‘correct’  choice  of  the  image  plane  position  has  been  the  subject  of  much 
debate.  It  is  not  deeper  into  the  metal  than  the  first  nuclear  plane  for  low 
index  faces.  This  is  the  position  we  will  assume.  For  a  grounded  ideal  metal, 
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the  relative  permittivity  tT  is  finite,  and  =  —q.  For  a  molecule  modeled 
by  point  sites,  each  point  charge  interacts  with  all  the  image  charges  by  means 
of  a  Coulomb  potential 

V„-„w  =  £  ^  (3) 

l,k  r* 

where  the  index  /  stands  for  the  real  charge  and  k  for  the  image  charge  with 
rjjt  being  their  distance  apart.  Since  the  polarization  response  of  the  surface 
changes  with  the  addition  of  new  molecules,  Eq.  (3)  describes  a  non-additive 
surface  polarization. 

The  second  term  in  Eq.  (1)  contains  the  anisotropic  short  range  potential 
and  the  short  range  r~n  potential  (n  integer)  and  so  requires  much  more 
detailed  description. 

core  ~  Vw— on  ~f"~  Lfi/—  ttr  (4) 

where  I4,_,v  is  the  isotropic  short  range  part  and  Vw~an  contains  the  anisotropic 
part.  Next  we  construct  both  these  parts  separately  using  pair-wise  sums  and 
then  use  translational  symmetry  to  reduce  them  to  usable  forms.  Henceforth 
we  drop  explicit  reference  to  water  (label  to) 

To  begin,  consider  a  perfect  crystal  surface  in  which  the  basic  vectors 
of  the  lattice  are  ai  and  a2.  Next  consider  a  point  particle  p  at  rp=(pp, 
zp)  where  zp  is  the  perpendicular  distance  above  the  surface  and  pp  is  the 
projection  of  rp  onto  the  surface  plane.  The  label  p  stands  for  0  or  either 
of  two  H  atoms.  For  convenience  we  choose  the  x-axis  to  be  parallel  to  aj 
and  take  origin  to  be  at  the  inversion  center  of  the  first  lattice  plane.  Let 
the  atoms  j  of  the  surface  plane  be  located  at 

1,  =  nxjaj  +  n2ja2  (5) 

then  the  distance  between  atom  j  and  the  particle  p  is  written 

r  .  —  -  / 


rw  =  y  P pi  +  zPi 

where  ppi  =  xpj  +  ypj-,  and  xpj,yw-,  zpj  are  the  components  of  the  vector 

rw  =  rp  —  h  (^) 

In  this  frame  the  anisotropic  part  of  the  potential  has  the  form 


v„tar„)  -  4vE  ((b  4)  4) 
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Here  a  is  a  scaling  parameter  that  acts  differently  on  the  repulsive  and  at¬ 
tractive  parts  and  gives  the  potential  its  anisotropy  between  directions  in  the 
xy  plane  and  the  z  axis.  In  our  calculations a  was  set  equal  to  0.8,  since  much 
smaller  values  would  not  be  consistent  with  weakly  chemisorbed  molecules. 
The  summation  in  Eqn.(8)  is  taken  over  all  atoms  j  of  the  surface. 

The  short-range  isotropic  contribution  can,  in  general,  be  described  as 

V..An,p-.r,)  =  (9) 

i  Pi 

where  Cn  is  a  fitting  parameter  with  n  being  an  integer.  The  scale  factors 
<Tpt  and  eps  are  separated  for  mathematical  convenience. 

Now  that  we  have  described  the  form  of  the  water-core  interactions  for 
the  two  models  we  can  write  down  the  complete  potentials  for  the  two  cases. 
In  the  first  model  (Al)  we  have 


Vax  =  Vw.cond  +  Van(0\  To)  +  Vl*r(8,0;  rG)  -  £  V{„{ 6,  H;  rH)  (10) 

H 


In  the  second  model  (A2)  we  have 


VA2  =  K,_w  +  Von(0;  to)  +  Vi„(  10, 0;  ro)  +  E[^n(H;  rH)  +  K„(10,  H;  rH)] 

H 

(if) 

Note  that  in  the  second  model  A2,  both  the  oxygen  and  hydrogen  atoms  in 
water  interact  with  the  metal  atoms  through  the  anisotropic  Lennard-Jones 
potential,  and  the  short-range  interaction  now  has  an  inverse  tenth  power 
dependence.  This  latter  change  is  needed  to  describe  the  atomic  polarization 
under  strong  electric  field  near  the  metal  surface.  The  physical  origin  of  such 
electrostatic  field  is  the  deformation  of  the  electron  cloud  within  a  metal 
[27,  28,  29].  We  shall  show  that  model  A2  improves  the  binding  energy,  the 
equilibrium  distance,  and  the  adsorption  via  hydrogen  atoms  over  model  Al 
but  at  the  cost  of  more  complexity. 

Following  the  derivation  of  Steele  [24,  25],  the  lattice  sum  in  Eqn.(8)  is 
transformed  to  reciprocal  space: 


^on(p;  rp) 


8*ep,ap2,  f  <7™  aV;t 

i4ceii  \  10a2  z10  4z4 
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£«<*•*  •  *•»  (I)1  *<**/•>  -  ^  (£)'**«.>] }  a*) 

In  the  above  equation,  the  label  p  is  omitted  on  the  right  hand  side  from  p 
and  z,  Aoeji  is  the  surface  area  of  a  unit  lattice  cell,  Km(£)  stands  for  mth 
order  modified  Bessel  function  of  the  second  kind,  and  G(gj,gj-  p)  includes 
contributions  from  all  the  trigonometric  factors  coe(gj-  p)  with  the  restriction 
of  gj  =  constant. 

Analogously  the  short  range  component  can  be  written 


Vif{n>Pirp)  — 


(t r + ?  (i)"'1" 


in  terms  of  the  same  set  of  geometric  factors  G.  Note  that  the  G  factors 
depend  on  the  symmetry  of  the  lattice  plane  and  are  different  for  (111), 
(100),  and  (110)  planes. 

We  conclude  with  explicit  lists  of  the  geometric  factors  for  each  of  the 
three  main  low  index  planes  starting  with  the  square  lattice.  For  the  (100) 
plane  where  |ai|  =  |a2|  —  afyj 2  and  A^u  =  a2/ 2,  the  lengths  of  the  two- 
dimensional  g j  vectors  in  the  reciprocal  lattice  space  are  given  by 


9>  =  ~  {3ji +9%)' 


where  a  is  the  lattice  constant,  g,\  and  g,2  are  a  set  of  integers,  which  satisfy 
gj  =  constant.  Accordingly,  we  write  the  dot  product  of  g j  and  p  as 


g  jP  =  — («#i*  +  9iiV)- 
a 


Straightforward  manipulation  of  Eqns.(14)  and  (15)  leads  to  the  following 
explicit  expressions  for  the  first  five  leading  terms  of  the  corrugation  factors 
G: 

G(gi;  gi  •  p)  =  cos  0xx  +  cos  0yy, 

G{g2\  g2  •  p)  =  2  cos  0xx  cos  0yy, 

Gigzl  g3  •  p)  =  cos  2  0xx  +  cos  20vy, 
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and 


G(S4‘>  g4  •  p)  =  2  (cos  2/9rx  cos  /9„y  +  cos  /9*x  cos  20vy ) , 


Gfas;  gs  •  p)  =  2  cos  4&x  cos  4/3„y, 

corresponding  to  gx  =  y2  =  %/2&,  y3  =  2/3r,  =  \/50x,  and  y5  =  2 y2. 

In  the  above  equations,  (3X  =  =  2\/2 x/a. 

The  situation  for  (110)  surface  is  similar,  except  for  |a2 1  =  |a2|/\/2  = 
a/y/2.  Accordingly,  equation  (14)  is  replaced  by 

&  =  ~(2&  +&Y'2  (16) 

which  leads  to 

gj  •  p  =  ~(y/2g,tx  +  gj2y). 

Since  the  Fourier  series  converges  relatively  slowly  at  small  z  for  this  surface, 
we  keep  the  first  ten  terms,  corresponding  to  g\  =  2 x/a,  g2  =  \/2gi,  g3  = 
y/%9i,  9*  =  2</i,  y5  =  S«  =  2y2,  y7  =  3y,,  g%  =  \/lTyi,  y»  =  2g3,  and 
<7io  =  Agi .  These  corrugation  components  have  the  form  of 

G($i;gi  p)  =  cos/3vy, 


G{g2,Zi  -  p)  =  cos/?rx, 
<7(33!  g3  •  p)  =  2  cos  (}xx  cos  /?yy , 
Gfoigs-p)  =cos2/3„y, 
C?(^5;  g5  •  p)  =  2  cos  /3xx  cos  2£„y, 
G(y6;gs  •  p)  =  cos2/3rx, 


8 


In  the  case  of  (111)  plane,  |ai|  =  |a2|  =  a/y/ 2  and  =  y/3a2/2.  The 
magnitude  of  gj  is  evaluated  by 

Si  =  (s*i  +  sh  ~  Susa)1'2 ,  (17) 

which  gives  51  =  2/9„,  y2  =  v/3$x,  y3  =  2yx,  54  =  \/7yi,  ffs  =  3yx  ••  •  with 
Px  —  2y/2r/a  and  @¥  =  2Px/y/3.  For  such  configuration,  the  dot  product 

2X 

gj  p=  —  [(2<7>i  - 9j2)x/V3  +  5ijy] ,  (18) 

yielding  the  first  five  terms  in  the  Fourier  series: 

G(gi ;  gi  •  p)  =  cos  2pyy  +  2  cos  pxx  cos  fiyy, 

G{g2]  g2  •  p)  =  cos  2/3xx  4-  2  cos  (3xx  cos  3 fivy, 

G(g3-,  g3  •  p)  =  cos  40„y  +  2  cos  2fixx  cos  2 /S„y, 

G(g< ;  g<  •  p)  =  2(cos  3ftTx  cos  f)vy  +  cos  2fixx  cos  40yy  +  cos  fixx  cos  5 fiyy) , 

G(g$]  g5  •  p)  =  cos  6 0yy  +  2  cos  3/3xx  cos  3/9yy. 

In  most  cases,  only  the  first  few  expansion  terms  are  necessary  to  achieve 
sufficient  accuracy.  However,  for  some  extreme  conditions,  there  may  be  a 
need  to  keep  more  terms. 

3  Results  and  Discussion 

In  this  section  we  describe  calculations  with  the  two  forms  of  the  potential 
for  different  low  index  planes  of  some  fee  metals.  In  all  the  calculations,  we 
set  a  =  0.8,  which  reduces  the  effective  radius  for  the  repulsion  and  increases 
the  radius  for  the  attraction  through  inverse  dependence  in  the  second  term 
of  Eqn.  (8),  when  the  projection  vector  pw  is  not  null.  The  simple  point 
charge  model  SPC/E  is  used  to  represent  the  water  molecule.  Its  geometry, 
partial  charges,  and  Lennard- Jones  parameters  are  used  for  model  Al.  For 
model  A2,  the  Lennard-Jones  parameters  were  replaced  by  taking  the  oxygen 
as  a  neon  atom  with  <ro  =  0.286nm,  Co  =  0.229kJ/mol  [30],  and  hydrogen  as 
a  helium,  with  <th  =  0.256nm,  and  ch  =  0.0845kJ/mol  [31]. 
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The  Lennard-Jones  parameters  for  the  metals  were  taken  from  the  lit¬ 
erature.  Table  I  lists  the  potential  parameters  evaluated  from  a  variety  of 
crystalline  state  physical  properties  at  a  wide-range  of  temperatures  [33,  34] 
and  the  lattice  constants  [30]  for  a  selected  set  of  metals.  To  get  the  coupling 
parameters  (so-called  mixed  Lennard-Jones  parameters)  needed  for  the  two 
forms  of  potential  described  in  the  last  section  we  used  the  parameters  just 
described  and  the  the  well-known  Kong’s  combining  rules  [32].  This  limits 
the  number  of  the  fitting  parameters  to  the  two  Cs  in  Eqns.(8  -  11),  except  in 
the  case  of  A1  where  there  was  only  one  parameter.  In  passing  we  comment 
that  since  the  Lennard-Jones  radii  and  energies  used  in  A2  do  not  depend 
on  the  particular  model,  readjustment  of  parameters  for  other  water  models 
that  have  different  geometry  and  point  charges  is  straightforward. 

Reported  in  Table  II  are  the  coefficients  for  the  angular  independent, 
short-range  interactions  in  Vi„  Eqns.(9  -  11),  obtained  by  a  fit  to  the  ab 
initio  binding  energies  of  water  on  Ni(lll)  [14],  Pt(100)  [10],  and  Al(100) 
[11].  For  the  case  of  aluminum,  we  simplified  the  calculation  further  by 
setting  Cio(H)  =  0.  In  this  study,  we  define  the  atop  site  as  the  position 
right  above  a  metal  atom,  the  bridge  site  as  the  position  above  the  middle 
point  of  the  ax  vector,  and  the  hollow  site  is  the  one  with  a  second  layer 
metal  atom  underneath.  There  is  a  complication  that  we  did  not  consider  in 
our  calculations.  According  to  the  calculations  of  Yang  and  Whitten  [14]  for 
Ni,  the  adsorption  energy  for  the  hollow  site  with  no  second  layer  metal  atom 
underneath  is  9.59  kJ/mol  larger  than  the  one  with  the  atom.  We  considered 
this  effect  to  be  beyond  the  scope  of  the  present  model  and  so  we  took  an 
average  value  of  the  two,  -28.65  kJ/mol,  as  the  criterion  for  comparison. 
It  should  also  be  mentioned  that  all  the  quantum-mechanical  calculation 
results  are  based  on  small  clusters,  and  may  depart  from  the  surface  results 
significantly. 

Although  the  fitting  procedure  is  for  specific  low  index  faces,  the  param¬ 
eters  determined  are  equally  applicable  to  other  types  of  crystallographic 
structure  as  well  as  to  discrete  atomic  layers.  This  can  clearly  be  seen  from 
Tables  III-V  where  comparisons  with  results  from  other  potential  functions 
and/or  molecular  orbital  calculations  are  given.  Keep  in  mind  that  two  as¬ 
sumptions  are  implicit  in  the  above  statement:  pair-wise  additivity  approx¬ 
imation  for  non-Coulombic  interactions  and  the  validity  of  the  combining 
rule.  Also  note  that  the  adsorption  energies  were  calculated  for  a  molecu¬ 
lar  orientation  with  the  oxygen  atom  closest  to  the  surface  (‘oxygen  atom 
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down’)  and  with  the  vector  bisecting  the  H-O-H  bond  angle  perpendicular 
to  the  interface,  and  the  line  connecting  the  two  hydrogen  atoms  parallel 
to  the  £*axis  (parallel  to  ai).  These  calculated  results  are  consistent  with 
experimental  measurements  {9,  35,  36,  37,  38]. 

According  to  our  calculations,  the  equilibrium  distances  between  the  oxy¬ 
gen  atom  of  the  water  molecule  and  the  Ni(lll)  surface  at  atop,  bridge,  and 
hollow  sites  are,  respectively,  0.249,  0.260,  0.261  nm  for  model  Al  and  0.226, 
0.240,  0.243  nm  for  model  A2,  compared  with  0.206,  0.209,  0.210  nm  of  Yang 
and  Whitten  [14].  Both  our  models  predict  an  equilibrium  distance  that  is 
too  large  compared  to  the  ab  initio  calculations.  Model  A2  is  marginally 
‘better’  than  Al  in  yielding  the  bond  distances  that  are  smaller  and  closer 
to  the  cluster  calculation. 

According  to  the  calculation  of  Yang  and  Whitten  [14],  the  adsorption 
energy  of  water  on  Ni(lll)  surface  is  shifted  upward  by  5.02  kJ/mol  when 
the  molecule  on  the  atop  site  rotates  about  the  C2V  axis  by  90°  so  the  H-H 
line  becomes  parallel  to  the  y-axis.  In  our  model,  the  energy  difference  turns 
out  to  be  much  smaller  since  the  distances  from  the  hydrogen  nuclei  to  the 
nearest  Ni  atom  remain  the  same  after  the  rotation.  It  seems  that  the  energy 
shift  observed  by  Yang  and  Whitten  arises  from  some  binding  contributions 
absent  in  the  present  model. 

It  is  informative  to  examine  the  energy  variation  with  respect  to  4>  which 
is  defined  as  the  angle  between  the  molecular  dipole  moment  vector  and  the  z 
axis  normal  to  the  surface.  In  particular,  rotating  the  water  molecule  on  the 
Ni(lll)  surface  (keeping  zq  fixed  at  the  equilibrium  distance)  by  <f>  =  49.75° 
in  the  xz  plane  increases  the  binding  energy  to  -4.19  kJ/mol  for  model  Al 
and  to  -9.52  kJ/mol  for  model  A2  at  atop  site,  to  -25.24  kJ/mol  and  18.51 
kJ/mol  at  the  bridge  site,  and  to  4.83  kJ/mol  and  -15.37  kJ/mol  at  the 
hollow  site.  On  the  other  hand,  the  value  calculated  by  Yang  and  Whitten 
[14]  for  the  atop  site  at  <f>  =  52.25°  is  25.10  kJ/mol.  There  is  a  5°  difference 
in  the  H-O-H  bond  angle  of  the  SPC/E  model  compared  to  the  geometry  of 
water  in  the  ab  initio  cluster.  However  we  have  chosen  <j>  such  that  one  of 
our  O-H  bonds  makes  the  same  angle  relative  to  the  z  axis  as  in  the  Yang- 
Whitten  calculation.  Clearly  our  relatively  simple  model  underestimates  the 
rotational  barrier. 

By  allowing  the  molecular  dipole  to  change  its  direction  along  the  yz 
plane,  using  xj)  to  measure  the  angle  between  the  dipole  vector  and  the  surface 
normal,  we  are  able  to  study  the  binding  energy  variation  upon  molecular  tilt. 
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Successively  tilting  the  water  molecule  in  increments  of  10°  from  the  surface 
normal  (ip  =  0)  changes  the  binding  energy  at  the  top  site  to  -43.82  (ip  =  10°), 
-42.20  (rP  =  20°),  -39.45  (V»  =  30°),  -35.51  (iP  =  40°),  -30.25  (ip  =  50°),  -23,57 
(ip  =  60°)  kJ/mol  for  model  A1  and  to  -44.32,  -44.02,  -43.50,  -42.77,  -41.76,  - 
40.32  kJ/mol  for  model  A2.  By  comparison,  the  corresponding  values  of  Yang 
and  Whitten  were:  -47.28,  -48.12,  -48.12,  -45.19,  -41.00,  -35.98  kJ/mol.  The 
energies  computed  from  the  simple  models  show  a  monotonic  increase,  the 
quantum  calculation  result  displays  a  minimum  around  ip  =  20°  to  30°.  These 
discrepancies  should  not  affect  the  computer  simulation  results  significantly 
since  they  only  affect  the  water  in  contact  with  the  surface  and  decrease 
rapidly  with  increasing  distance  z. 

The  Lennard-Jones  parameters  obtained  from  thermodynamic  data  sug¬ 
gest  that  Cu,  Ag,  Au,  and  Ni  form  a  single  class,  and  Pt  yet  smother  with 
stronger  binding.  Using  the  fitting  coefficients  of  Ni,  we  are  able  to  predict 
the  energies  for  adsorption  on  copper  surfaces.  The  results  axe  given  in  Ta¬ 
ble  VI,  and  show  that  the  predicted  binding  energy  of  water  molecule  right 
above  a  metal  atom  of  Cu(100),  -38.70  kJ/mol  at  0.251  nm  for  model  A1 
and  -38.95  kJ/mol  at  0.228  nm  for  model  A2,  is  in  accord  with  the  ab  initio 
results  of  Ribarsky  and  co-workers  [12},  namely  -36.66  kJ/mol  at  0.220  nm. 

All  the  potential  surfaces  depicted  in  Tables  III-VI  indicate  a  preference 
for  water  to  adopt  orientations  in  which  water  dipole  moment  points  away 
from  the  interface  and  the  oxygen  sits  down  on  the  surface  with  the  on  top 
site  being  the  site  with  greatest  binding  energy.  If  more  experimental  data 
were  available  the  potentials  described  here  could  be  further  refined. 

Using  the  parameters  reported  in  Tables  I  and  II,  we  plot  the  potential 
energy  curves  calculated  from  the  analytical  formulas  derived  in  Section  2 
for  one  oriented  water  molecule  adsorbed  at  different  sites  on  low  index  faces 
of  platinum.  The  results  are  shown  by  the  curves  in  Figures  1-3.  For  com¬ 
parison,  we  also  plot  the  potential  energies  (data  points)  obtained  from  the 
direct  summation  of  pair-wise  interaction.  Both  sets  of  calculations  explic¬ 
itly  included  the  top  two  layers  of  metal  surface  atoms  and  approximated 
the  remainder  of  the  crystal  by  means  of  the  9-3  potential  that  comes  from 
treating  the  remaining  half  spaue  in  the  continuum  limit.  As  can  be  seen 
there,  the  agreement  between  these  two  sets  of  computations  is  excellent. 
The  largest  deviations  hardly  visible  in  the  plots  occur  when  the  particle- 
wadl  distance  is  very  short.  Since  the  use  of  the  transformed  form  of  the 
potential  energy  is  computationally  much  more  economical,  direct  summa- 
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tion  of  pair-wise  interaction  is  not  necessary  in  computer  simulations  unless 
inclusion  of  the  thermal  motions  of  the  metal  atoms  is  desired.  In  Figs.  1-3, 
about  75%  of  the  adsorption  energy  on  atop  site  comes  from  the  interac¬ 
tion  with  the  nearest  metal  atom.  For  bridge  site,  on  the  other  hand,  the 
contribution  from  the  two  nearest  neighbors  is  about  65%.  By  comparing 
the  potentials  at  distances  greater  than  0.35  nm  in  Figures  1-3  we  see  that 
the  corrugation  potential  decays  very  rapidly  with  distance  z  from  the  metal 
surface.  However  if  the  first  layer  of  water  is  strongly  bound  in  an  ordered 
array  then  effects  of  corrugation  will  be  felt  further  out. 

When  angle  dependence  of  the  potential  is  turned  off  (or  =  1  and  C„(p)  = 
0),  the  top  site  becomes  the  site  of  weakest  binding.  This  is  shown  for  the 
Pt(lll)  surface  and  model  A2  in  Figure  4.  Figure  4  shows  how  anisotropy 
is  important  in  lowering  energy  of  the  top  site  and  in  decreasing  the  metal- 
oxygen  bond  distance.  As  mentioned  above,  introduction  of  the  angular 
dependent,  short-range  interactions  to  represent  chemical  bonding  is  essential 
for  mimicking  water  adsorption  on  metallic  surfaces.  The  Lennard-Jones 
potential  function  does  not  produce  correct  preferential  adsorption  site.  By 
adding  the  angular  dependent,  short-range  interactions,  the  potential  well  for 
atop  site  deepens  and  shifts  towards  small  z,  while  for  the  bridge  and  hollow 
sites  the  influence  is  just  opposite.  This  reverses  the  preference  order  for 
adsorption  energies.  The  contribution  from  image  interactions  is  also  shown 
in  Figure  4. 

In  Figure  5  we  plot  the  total  energy  of  H2O-Pt(100)  system  versus  the 
O-surface  separation.  The  molecular  orientation  is  such  that  one  O-H  bond 
vector  is  perpendicular  to  the  horizontal  surface  plane  with  one  hydrogen 
pointing  directly  at  the  surface.  The  other  O-H  bond  lies  in  the  xz  plane 
pointing  away  from  the  surface  at  an  oblique  angle.  This  mimics  the  situation 
where  chemisorption  is  through  the  hydrogen  atom.  For  such  a  configuration, 
we  obtain  the  binding  energies  of  -26.11  kJ/mol  at  0.334  nm  for  the  atop  site, 
-15.82  kJ/mol  at  0.349  nm  for  the  bridge  site,  and  -13.38  kJ/mol  at  0.352  nm 
for  the  hollow  site  (model  A2).  The  corresponding  values  from  the  molecular 
orbital  calculation  of  Holloway  and  Bennemann  [10]  are  -24.18  kJ/mol  at 
0.265  nm  for  the  first  position  and  -5.94  kJ/mol  at  0.268  nm  for  the  third 
position.  Model  A1  does  not  predict  any  significant  adsorption  via  hydrogen 
atom.  Note  that  the  shift  in  minima  for  oxygen  versus  hydrogen  down  is 
approximately  0.1  nm.  This  is  the  shift  observed  by  Toney  et  al  [4]  in  the 
experiments  with  bulk  water  next  to  a  silver  surface  when  the  electric  field 
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is  reversed.  This  is  an  added  reason  for  adopting  potentials  with  the  form 
described  here. 

Compared  with  other  analytical  potential  functions,  the  present  one  is 
physically  more  meaningful.  It  contains  a  small  enough  number  of  fitting  pa¬ 
rameters  to  not  overly  prejudice  the  physics.  In  addition,  the  functional  form 
is  universal  for  studying  water  adsorption  on  various  metal  surfaces  and  can 
be  directly  used  for  systems  in  which  the  metal  atoms  are  explicitly  taken 
into  account.  As  new  experimental  or  theoretical  data  on  the  adsorption 
energy,  etc.,  becomes  available,  it  will  be  very  easy  to  modify  the  potential 
function  and  refine  the  parameters.  One  possible  channel  of  improvement  of 
the  present  model  is  to  add  non-additive  interactions  such  as  electric  polar¬ 
izations  or  other  many-body  interactions  explicitly.  Work  on  modifying  the 
potentials  to  include  electric  polarization  is  in  progress. 
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Table  I  Potential  parameters  [33]  and  lattice  constants  [30]  for 
metals. 


Metals 

e  (kJ/mol) 

a  (nm) 

a  (nm) 

A1 

37.84 

0.262 

0.405 

Cu 

39.49 

0.234 

0.362 

Ni 

mmm 

0.228 

0.352 

Pt 

65.77 

0.254 

0.392 

various 
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•  *♦ 


Table  II  Potential  parameters. 


Metals 

Cio  (0) 

Cw  (H) 

Cs(  0) 

Ce(H) 

Al 

1.56 

0 

2.16 

4.4 

Ni 

1.195 

0.8 

1.603 

3 

Pt 

1.25 

1.2 

1.331 

1.35 

a  =  0.8,  ow  =  0.316  nm,  oo  =  0.286  nm,  oh  =  0.256  nm, 
kJ/mol,  to  =  0.229  kJ/mol,  ch  =  0.0845  kJ/mol. 


=  0.650 
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Table  III  Adsorption  energies  of  water  on  nickel  surfaces  (in  units  of 
kJ/mol). 


A1 

A2 

YW 

SS 

- 

-49.0 

(100)  B 

-27.68 

-28.56 

- 

- 

(100)  H 

-22.08 

-23.65 

- 

- 

(110)  A 

-43.12 

-42.85 

- 

-55.9 

(110)  B 

-40.46 

- 

- 

(110)  H 

-26.77 

- 

- 

(111)  A 

-44.36 

-44.43 

-44.35 

-43.5 

(HI)  B 

-29.74 

-29.85 

-29.71 

-28.3 

(HI)  H 

-28.69 

-28.28 

-28.65 

-27.4 

A:  A-top;  B:  Bridge  site;  H:  Hollow  site;  ZP:  Present  work;  YW:  Ref. 
[14];  SS:  Ref.  [19];  Same  for  Tables  IV- VI. 


Table  IV  Adsorption  energies  of  water  on  aluminum  surfaces  (in  units  of 
kJ/mol). 


A1 

A2 

MH 

rtioiw.ii 

-51.06 

-51.27 

-51.14 

(100)  B 

-22.69 

-23.12 

-24.12 

(100)  H 

-15.17 

-17.86 

-14.47 

(110)  A 

-51.32 

-50.58 

- 

(110)  B 

-22.80 

-22.01 

- 

(110)  H 

-8.36 

-13.89 

- 

(111)  A 

-50.78 

-51.27 

- 

(HI)  B 

-24.29 

-24.08 

- 

(HI)  H 

-22.83 

-22.39 

- 

MH:  Ref.  [11] 
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Table  V  Adsorption  energies  of  water  on  platinum  surfaces  (in  units  of 
kj/mol). 


A1 

A2 

SS 

SH 

RFMB 

HB 

FRB 

I B 

-48.25 

-48.23 

EEftl 

- 

-48.25 

-39.5 

(100)  B 

-33.18 

-28.88 

-22.2 

-16.3 

- 

- 

-22.0 

(100)  H 

-27.94 

-20.7 

-12.1 

- 

-27.94 

-8.2 

(110)  A 

-46.03 

-47.28 

- 

- 

- 

- 

- 

(110)  B 

-30.96 

-27.39 

- 

- 

- 

- 

- 

(110)  H 

-20.72 

- 

- 

- 

- 

- 

(111)  A 

-49.74 

Ena 

- 

- 

-40.7 

- 

- 

(HI)  B 

-35.67 

-30.17 

- 

- 

-26.1 

- 

- 

(HI)  H 

-34.57 

-28.12 

- 

- 

-23.9 

- 

- 

SH:  Ref.  (16];  RFMB:  Ref.  (18];  HB:  Ref.  (10];  FRB:  Ref.  [17]. 
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Table  VI  Adsorption  energies  of  water  on  copper  surfaces  (in  units  of 
kJ/mol). 


A1 

A2 

RLL 

Ullil 

-38.70 

imm 

- 

(100)  H 

-18.97 

- 

(110)  A 

-38.20 

-38.16 

- 

(110)  B 

-23.25 

-23.85 

- 

(110)  H 

-12.91 

-16.78 

- 

(111)  A 

-39.12 

-39.52 

- 

(HI)  B 

-26.19 

- 

(HI)  B 

-24.73 

-24.76 

- 

RLL:  Ref.  [12] 
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Figure  Caption 

Fig.  1  Water-Pt(lOO)  potential  (model  A2)  as  a  function  of  distance  of 
oxygen  atom  from  surface  plane.  The  molecule  lies  in  y  =  constant 
plane  with  dipole  moment  vector  parallel  to  surface  normal  and  oxygen 
atom  down.  Symbols:  direct  sums  of  pair-wise  energy.  Curves:  Fourier 
series  expansion. 

Fig.  2  Water-Pt(IlO)  potential  (model  A2)  as  a  function  of  distance  of 
oxygen  atom  from  surface  plane.  The  molecule  lies  in  y  =  constant 
plane  with  dipole  moment  vector  parallel  to  surface  normal  and  oxygen 
atom  down.  Symbols:  direct  sums  of  pair-wise  energy.  Curves:  Fourier 
series  expansion. 

Fig.  3  Water-Pt(lll)  potential  (model  A2)  as  a  function  of  distance  of 
oxygen  atom  from  surface  plane.  The  molecule  lies  in  y  =  constant 
plane  with  dipole  moment  vector  parallel  to  surface  normal  and  oxygen 
atom  down.  Symbols:  direct  sums  of  pair-wise  energy.  Curves:  Fourier 
series  expansion. 

Fig.  4  Water-Pt(lll)  potential  without  anisotropic  and  short-range  isotropic 
terms.  The  potential  is  a  sum  of  the  conventional  12-6  Lennard- 
Jones  interaction  and  the  Coulomb  interaction  between  real  and  image 
charges.  The  Coulombic  image  potential  is  also  plotted  separately  to 
show  its  weak  z  dependence.  The  molecule  lies  in  y  =  constant  plane 
with  dipole  moment  vector  parallel  to  surface  normal  and  oxygen  atom 
down. 

Fig.  5  Water-Pt(lOO)  potential  (model  A2)  with  H  atom  down.  The  poten¬ 
tial  is  plotted  as  a  function  of  distance  of  oxygen  atom  from  the  surface 
plane.  The  molecule  lies  in  y  —  constant  plane  with  one  0-H  bond  per¬ 
pendicular  to  the  surface.  Note  shift  in  minima  by  approximately  0.1 
nm  and  reduced  binding  energies  compared  to  the  oxygen  atom  down 
orientation  shown  in  Figure  1. 
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Fig.  1  Water-Pt(IOO)  interaction 
Symbols:  direct  sums  of  pair-wise  energy 
Curves:  Fourier  series  expansion  of  energy 
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Fig.  2  Water-Pt(IIO)  interaction 
Symbols:  direct  sums  of  pair-wise  energy 
Curves:  Fourier  series  expansion  of  energy 
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Fig  3  Water-Pt(lll)  interaction 
Symbols:  direct  sums  of  pair-wise  energy 
Curves:  Fourier  series  expansion  of  energy 


Potential  Ener 


Distance  from  Surface  (nm) 

Fig.  5  Water-Pt{100)  interaction 
Adsorbed  via  hydrogen 


